#include "matrix_tool.h"

arma::mat toeplitz_inverse(arma::mat mat)
{
	int dim = mat.n_cols;
	arma::vec a(dim-1);
	arma::vec b(dim-1); b.zeros();
	
	for (int i = 0; i < dim-1; ++i)
	{
		a(i) = mat(0,i+1)/mat(0,0);
	}
	// a.print("a=");
	b(0) = -a(0);
	for (int i = 1; i < dim-1; ++i)
	{
		for (int j = i-1; j >= 0 ; --j)
		{
			b(i) -= a(i-1-j)*b(j);
		}

		b(i) -= a(i);
	}
	b.insert_rows(0,arma::eye(1,1));
	// b.print("b=");
	b = b/mat(0,0);
	// b.print("b'=");
	return toeplitz(b,arma::zeros<arma::vec>(dim)).t();
}

arma::vec toeplitz_inverse(arma::vec vec)
{
	int dim = vec.n_rows;
	arma::vec a(dim-1);
	arma::vec b(dim-1); b.zeros();
	
	for (int i = 0; i < dim-1; ++i)
	{
		a(i) = vec(i+1)/vec(0);
	}

	b(0) = -a(0);
	for (int i = 1; i < dim-1; ++i)
	{
		for (int j = i-1; j >= 0 ; --j)
		{
			b(i) -= a(i-1-j)*b(j);
		}

		b(i) -= a(i);
	}
	b.insert_rows(0,arma::eye(1,1));
	b = b/vec(0);
	return b;
}

arma::vec col_toeplitz_mul(arma::vec toeplitz_vec, arma::vec vec)
{
	arma::vec x(vec.n_rows);
	x.zeros();
	for (int i = 0; i < vec.n_rows; ++i)
	{
		for (int j = i; j >= 0; --j)
		{
			x(i) += toeplitz_vec(j)*vec(i-j);
		}
	}

	return x;
}

arma::vec row_toeplitz_mul(arma::rowvec toeplitz_row, arma::vec vec)
{
	arma::vec x(vec.n_rows);
	x.zeros();
	for (int i = 0; i < vec.n_rows; ++i)
	{
		for (int j = 0; j <= vec.n_rows-i-1; ++j)
		{
			x(i) += toeplitz_row(j)*vec(j+i);
		}
	}

	return x;
}